%% plot CL CD
alphadeg = 10:0.1:20;
Ma = 10;
%% lidongxu
% CL = -0.21 + 0.005059*Ma + 0.0497*alphadeg + 0.0004244*alphadeg.^2;
% CD = 0.05244 + 0.00209*Ma - 0.01629*alphadeg + 0.001937*alphadeg.^2;
%% gyh
% p00 = 0.5956;  p10 = -0.09572; p01 = -0.06861; p20 = 0.0005887;
% p11 = 0.01532; p02 = 0.003406; p30 = 0.0005817; p21 = -0.00101;
% p12 = -0.000579; p40 = -3.246e-05; p31 = 2.228e-05; p22 = 3.615e-05;
% p50 = 5.154e-07; p41 = -5.753e-08; p32 = -7.175e-07;
% x = Ma;
% y = alphadeg;
% CD = p00 + p10*x + p01*y + p20.*x.*x + p11.*x.*y + p02.*y.*y + p30.*x.*x.*x + p21.*x.*x.*y...
%      + p12.*x.*y.*y + p40.*x.*x.*x.*x + p31.*x.*x.*x.*y + p22.*x.*x.*y.*y + p50.*x.*x.*x.*x.*x...
%      + p41.*x.*x.*x.*x.*y + p32.*x.*x.*x.*y.*y;
% p00 = -0.04928;		p10 = -0.009818;	p01 = 0.05259;		p20 = 0.0011;
% p11 = -0.001381;	p02 = 0.0004019;	p30 = -2.586e-05;	p21 = 2.176e-05;
% p12 = -2.289e-06;
% CL = p00 + p10.*x + p01.*y + p20.*x.*x + p11.*x.*y + p02.*y.*y + p30.*x.*x.*x +...
%     p21.*x.*x.*y + p12.*x.*y.*y;
%% rlv
% CL = -0.041065 + 0.016292*x + 0.00026024.*x.*x;
% CD = 0.080505 - 0.03026*CL + 0.86495*CL.^2;
%% 乔浩
p1 = -0.04701; p2 = -0.006756; p3 = 0.05069; p4 = 0.0004161; p5 = -0.0008746; p6 = 0.0003743; 
q1 = 0.1672; q2 = -0.01845; q3 = 0.001565; q4 = 0.0007091; q5 = -0.0003823; q6 = 0.00083;
CL = p1 + p2 .* Ma + p3 .* alphadeg + p4 .* Ma.^2 + p5 .* alphadeg .* Ma + p6 .* alphadeg.^2;
CD = q1 + q2 .* Ma + q3 .* alphadeg + q4 .* Ma.^2 + q5 .* alphadeg .* Ma + q6 .* alphadeg.^2;
%% Terminal Guidance for a Hypersonic Vehicle with Impact Time Control
CL = 0.4172 + 19.41*alpha + 10.17*alpha^2 - Ma*(0.1004 + 0.7536*alpha);
CD = 0.3042 + 0.02988*CL^2; 
%% figure
figure
plot(alphadeg, CL, 'linewidth', 1.5);
xlabel('alpha'); ylabel('CL');
figure
plot(alphadeg, CD, 'linewidth', 1.5);
xlabel('alpha'); ylabel('CD');